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1. Introduction 

Present knowledge on the large scale infrared emission from 
interstellar dust is based largely on the all-sky surveys per- 
formed by the IRAS and the COBE satellites. IRAS was the 
first to observe the sky at a resolution of a few arc minutes at 
wavelengths between 12 and 100/im. Clear variations were ob- 
served in the infrared spectrum of interstellar clouds, e.g. the 
ratio of 60/im vs. 100/im dropping towards dense clouds (e.g. 
Laureijs et al. 1996| ; Abergel et al. 1994). These were inter- 
preted as the result of abundance variations between classical, 
large grains and very small dust particles that are transiently 
heated above the equilibrium temperature (e.g. Desert et al. 



199C). Emission in the mid-infrared and especially the IRAS 



12/im band, are now believed to be caused by even smaller 
grains or very large molecules, the most common candidates 
being Polycyclic Aromatic Hydrocar bons (PAHs) (Leger et al. 
|l989f ; Desert et al. |l990t Li & Draine |200l| ). The PAH emission 
traces the warmer outer surfaces of the clouds. 

The COBE/DIRBE instrument covered a wider wavelength 
range from near-infrared up to 240 /im and the FIRAS instru- 
ment up to 1 cm. The spatial resolution was, however, com- 
paratively poor i.e. -40' for DIRBE and -7° for FIRAS. 
Nevertheless, the extended wavelength coverage revealed sur- 
prisingly large quantities of cold dust that were not observed 
with IRAS. Detailed analysis of the COBE data has shown 
that at least two emission components are needed to explain 
the observed far-infrared spectra. For most of the diffuse re- 
gions the dust emission follows a modified b lack b ody law with 
temp erature T — 17.5 K (B oulanger et al. 1996t Dwek et al. 
1997 ). Lagache et al. ( 1998| ) found, however, a colder compo- 
nent (T — 15 K) associated with molecular clouds. Finkbeiner 



et al. ( 1999 ) reached similar conclusions but included an even 
colder, T ^9 K, component. 



Recent observations with the PRONAOS balloon borne 
observatory have confirmed the lowering of the colour tem- 
perature towards molecular clouds and, e.g. in the case of 
the Polaris flare {Ay —1), towards moderately dense re- 



gions (Stepnik et al. [200 1| ). Similar effects have been seen 
in ISOCAM studies towards cirrus like clouds (Ay -^0.5; 
Miville-Deschenes et al. 2002 ) where the temperature varia- 
tions caused by extinction are quite insignificant. The observed 
increase in the sub-mm emission can be explained by grain 
growth and the formation of large dust aggregates (Stepnik 
et al. 2001; Cambresy et al. 2001) and in dense regions the 
very small grains seem to have disappeared almost completely. 
Conversely, in the more diffuse medium, large quantities of 
small grains can be produced by grain shattering (see e.g. 
Miville-Deschenes et al. |2002|). 



Based on these observations and theoretical studies we 
have obtained a rough picture of the dust properties and in- 
frared emission from interstellar clouds. In the outer layers we 
have the PAH emission and the warm, very small dust grains. 
Deeper in the clouds the extinction of short wavelength ra- 
diation reduces mid-infrared emission. At the same time, the 
grains start to form larger aggregates and/or acquire mantles 
as gas molecules freeze onto them. This increases the far- 
infrared and sub-mm emissivity relative to the absorption at 
short wavelengths and the physical grain temperatures also 
drop. The colour temperature therefore drops to approximately 
10 K (or possibly even lower) in dense cloud cores (e.g. Juvela 
et al. |2002|). 



Theoretical and laboratory studies have provided us with 
some information on the properties of likely dust grain mate- 
rials. Together with radiative transfer calculations, these have 
made it possible to determine size distributions for the differ- 
ent dust components consistent with the emission properties 



(e.g. Desert et al. 1990; Li & Draine 2001). On flie oflier hand 
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radiative transfer models can be used to study the variations 
of these properties in individual objects. These previous exam- 
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pies show that the dust emission spectrum is variable over the 
whole range from near-infrared to sub-mm. The changes are 
related to the local radiation field and the local density. So far, 
modelling has concentrated either on the determination of the 
dust emission under given radiation conditions, or at most in 
one-dimensional, spherical models divided into some tens of 
cells. 

In the related field of the modelling of molecular line emis- 
sion, three-dimensional models have already been used for 



some time (P ark & Hong |1995| ; Park e t al. |1996| ; Juvela |1998 
Padoan et al. 119981 Juvela et al. 120011). In that field, inhomo- 



geneous source structure has an even stronger effect on the 
emerging radiation and, on the other hand, direct calculations 
are already possible for models with up to ^^100"^ cells or, with 



approximate methods, even higher (Ossenkopf 2002) 



The continuum emission is not as directly linked with the 
local, physical conditions. It is, nevertheless, affected by the 
variations in the radiation field (extinction of the external radia- 
tion field and internal sources) and the changes in the dust prop- 
erties. So far dust emission has not been studied with general 
three-dimensional cloud models. The main reason is the com- 
plexity caused by the transiently heated small dust particles. In 
order to be able to predict the emission accurately, the distribu- 
tion of the grain temperatures must be first determined in each 
computational cell. Each dust population must be discretized 
into separate grain size intervals. For each of these, a couple 
of hundred enthalpy bins are needed and the solution of the as- 
sociated set of linear equations takes typically several seconds. 
This becomes significant when the model contains millions of 
cells and, in practice, it limits direct calculations to one- or pos- 
sibly two-dimensional models. 

We will discuss how this limitation can be overcome by 
suitable approximations without siginificantly affecting the ac- 
curacy of the results. This will be based on the assumption that 
the local spectral energy distribution can be deduced from the 
intensity at suitably selected reference wavelengths. One must 
first determine a mapping between these intensities and the re- 
sulting local spectral energy distribution. The radiative transfer 
problem needs only to be solved at the reference wavelengths, 
and by using the established mapping, one can solve the radia- 
tive transfer problem. The solution is approximative but can be 
be made sufficiently accurate within the limitations of present 
day computers. The method makes it possible to study explic- 
itly the effects that the cloud structure has on the infrared emis- 
sion, and the method could be useful e.g. in the simulation of 
the spectral energy distribution of galaxies. 

We first explain the implementation of the radiative transfer 
program (Sect. ^ and discuss its accuracy (Sect. |]). In Sect. ^ 
we apply the method to three-dimensional model clouds that 
are based on numerical simulations of supersonic magneto- 
hydrodynamic (MHD) turbulence and we show some observ- 
able consequences of the inhomogeneity of the clouds. Finally, 



2. The program 

2. 1. Direct solution 

The full radiative transfer problem can be divided into two 
largely independent parts. In the first one determines the in- 
tensity of the radiation field at each position of the cloud (in 
our case in each discrete cell of finite volume). In the second 
part this information and the selected dust model are used to 
derive the local dust temperature and emission. In interstellar 
clouds dust emission mainly takes place in the infrared part of 
the spectrum. As long as the external radiation field is relatively 
weak the dust emission has only a small effect on the heating 
of grains (see e.g. Hollenbach et al. 1991 ; Bernard et al. 1992 ). 
This is true for all models considered in this paper and, there- 
fore, the two parts of the computation mentioned above can be 
carried out sequentially without iteration. In the following we 
will ignore other processes affecting the dust temperature. For 
example, the gas-dust coupling becomes important only at den- 



sities above JjlO cm (e.g. Kxiigel & Walmsley 1984) and 



in Sect. 4.3 we discuss some qualitative effects of spatial vari- 



ations in the dust properties. 



the clouds considered here have much lower densities. 

Our three-dimensional model clouds are divided into cubic 
cells, each with a constant density. First, the optical depths for 
absorption and scattering as well as the parameters of the scat- 
tering function are calculated for each cell. The usual Henyey- 
Greenstein formula (1941) and the asymmetry parameter g are 
used to define the scattering function. The radiation field is cal- 
culated with Monte Carlo methods. Model photons are initiated 
at random locations at the cloud boundary in order to simu- 
late the background radiation, and within the cloud to simu- 
late emission from the dust and possible internal sources. The 
path of the model photon is tracked in the cloud. A model pho- 
ton represents a number of real photons. Absorbed photons are 
counted in each cell that is crossed, and from time to time, 
the model photon is scattered towards a new direction as deter- 
mined by the dust model. Because of the frequency dependence 
of the scattering probability, the simulation has to be carried out 
separately for each frequency. 

The simulation determines the spectrum of the incoming 
radiation for each cell. In the second phase this is used to com- 
pute the dust temperature. The dust model consists of several 
dust populations, and each population is further discretized ac- 
cording to grain size. For each of these we have a probability 
distribution of the grain enthalpy, discretized usually to a few 
hundred bins. The transition probabilities for cooling processes 
are calculated beforehand, based on the specific heat and ab- 
sorption cross sections. The transitions upwards in energy are 
obtained by integration from the simulated intensity spectrum. 
The formulae we use are exact only for a probability distri- 
bution that is flat a function of grain energy (temperature). We 
have also implemented the 'thermal discrete' method described 
by Draine & Li (2001 ) and have tested that our method gives, 
with the used discretizations, identical results for the studied 
wavelength range. Once the number of grains as the function 
of temperature (enthalpy) is known the local emission can be 
calculated. 

The two steps still need to be iterated, if the dust emission 
itself contributes to the grain heating. This requires, however. 
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either the presence of strong, local heating sources or very high 
extinction. Bernard et al. (1992) studied the case of spherical 
clouds heated by normal ISRF (a slightly modified version of 
the ISRF given by Mathis et al. ( |1983[ )). The IR heating was 
found to be insignificant for visual extinctions Ay < 100™ 
through the cloud, and even with Ay — 1000'", it increased the 
sub-mm emission only by ^5%. The effect is due to an increase 
in the temperature of the large, classical grains. Shorter wave- 
length emission is in practice not affected, since the infrared 
photons emitted by the dust grains are not energetic enough to 
cause significant transient heating of the smaller grains. 

2.2. The library method 

For large, three-dimensional models consisting of > 100'^ cells 
the described direct approach is no longer possible. Simulation 
of the necessary number of frequencies (> 100) could still be 
done, but it would possibly take days on a workstation, depend- 
ing on the required accuracy. The main problem is the solution 
of the linear equations. Each dust population is discretized into 
ten or more size intervals, each further divided into possibly 
several hundred enthalpy intervals. The solution of the asso- 
ciated linear equations for the transition probabilities requires 
more than a second per cell, or a month for a model with 100'^ 
computational cells. The direct approach is clearly limited to 
very small models (or very large supercomputers). 

The spectrum of the local radiation field is very similar in 
many parts of the cloud. For example, in a homogeneous, semi- 
infinite plane or a homogeneous sphere, the intensity and spec- 
trum of incoming radiation depends only on one parameter, the 
distance from the free surface. In these cases one could calcu- 
late the dust emission at certain depths and obtain the emission 
from intermediate positions with interpolation. 

Inside an inhomogeneous or an iiTegularly shaped cloud 
there is more variation in the spectral energy distribution of the 
incoming radiation. Usually, in a homogeneous cloud, the ex- 
tinction reduces radiation more strongly at short wavelengths. 
This is not true, however, in the extreme case of inpenetrable 
clumps. If clumps absorb all incoming radiation irrespective of 
its wavelength, the intensity inside the cloud will be reduced 
proportionally to the solid angle blocked by the clumps, while 
the shape of the spectrum will remain unchanged. Clearly, in 
a general case, the changes in the radiation field cannot be 
tracked with just one parameter (e.g. visual extinction). We ar- 
gue, however, that for a given model, the spectrum can be suffi- 
ciently constrained with knowledge of the intensity at possibly 
only two or three reference wavelengths, Xf^^. 

The first step in our approximation is to find the range of 
intensities and spectral energy distributions of the incoming ra- 
diation. For this purpose, we reduce the discretization of the 
model down to ~ 32^ cells. This allows fast solution of the 
radiative transfer part of the problem. The result is the sim- 
ulated intensity in each cell. We take the range of intensities 
found at the first reference wavelength, A^°^, and discretize that 
into logarithmic intensity intervals. The procedure is repeated 
recursively at the other reference wavelengths. The result is a 
tree structure that represents a discretization of all observed in- 





Fig. 1. A schematic view of the tree structure used in the li- 
brary method. At each level the intensity of the incoming radi- 
ation at one wavelength is examined and each subnode cor- 
responds to a different intensity interval, Al. The leaves of 
the tree correspond to a unique combination of intensity values 
at the reference wavelengths and they hold the corresponding 
emission spectrum. Once the intensity of the incoming radia- 
tion is known at reference wavelengths, the emission can be 
read from the tree. 



tensity combinations, (/(A™^), I(A2°^), . . . ). In the leaves of the 
tree we store the average, full spectrum that corresponds to a 
unique combination of intensities at the reference wavelengths. 
If the reference wavelengths are suitably selected, the tree and 
knowledge of intensities at reference wavelengths can be used 
to obtain an estimate of the full spectrum. 

The second step of our approximation is the solution of the 
dust emission corresponding to each spectrum of incoming ra- 
diation that was stored in the leaves of the tree. This way, we 
obtain a mapping from the incoming intensity at the reference 
wavelengths to the coiTesponding full spectrum of local dust 
emission. We call this structure 'the library' (see Fig. |l]). 

To predict dust emission from a large, three-dimensional 
model, we first simulate the radiation field at the small number 
of reference wavelengths, Xf^ . The radiative transfer equation 
can then be directly integrated through the cloud, using in each 
cell the library to transform the simulated values of incoming 
intensity at the reference wavelengths into the full spectrum of 
emitted intensity. We call this 'the library method'. When we 
look up the emission spectrum from the library, we also per- 
form linear interpolation based on the values /(A[''^). Again, 
if absorbed dust emission is important (e.g. close to strong 
sources or in very opaque regions) iteration is needed. 
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3. Tests of the library method 

3. 1. The model clouds and the dust model 

We will now estimate the accuracy of the approximate method 
described in Sect. The method is based on the discretiza- 
tion of the incoming intensity into bins for which the cor- 
responding dust emission spectrum is only computed for the 
average spectrum. The accuracy can therefore always be im- 
proved by selecting a finer discretization and/or including more 
reference wavelengths. In order to obtain a significant advan- 
tage over the direct solution, the total number of spectrum bins 
must, however, be much less than the number of cells in the 
model. 

In this paper, we limit the study to starless clouds heated by 
the interstellar radiation field. We take as a test case an MHD 



simulation which is the model C from Juvela et al. (2001 ). The 



model is both supersonic and super- Alfvenic with acoustic and 
Alfvenic Mach numbers M « Ma wlO. The MHD simulation 
included self-gravity. 



The dust model is taken from Li & Draine (2001 ) and in 



eludes silicate grains (sizes a>3.5A), graphite grains (a ^50A) 
and PAHs (from a=3.5Ato ^50A). The PAH properties corre- 
spond to cold neutral medium. The grain size distributions are 



the same as in Li & Draine (2001). The only difference is the 



way PAH and graphite grain populations are joined. Instead 
of a gradual change in the description of absorption coefficient 



(Li & Draine 2001 , Eq. (2)) we only join the grain size distribu- 



tions smoothly around 40/im. The optical properties are taken 
from the files available through the website of B. Draine. The 
original references are Laor & Draine ( 1993 ) for graphite and 
silicate grains and Li & Draine (2001) for the PAH properties. 

At this point, we are interested in the differences between 
the two methods (i.e. direct solution and the library method) 
and not in the absolute accuracy. Therefore, the number of 
grain size intervals (~10), grain energies (<300) and simulated 
frequencies (^100) wer e kep t to a minimum that still repro- 
duces the Li & Draine ( 2001 ) spectra for diffuse ISM within 
~10% accuracy over the studied wavelength range, 5-300 /im. 
The largest differences result from a different method used to 
join the PAH and graphite grain populations around grain size 
40/ im. Th e interstellar radiation field adopted is from Mathis et 
al. (|1983|). 



The original cloud model represented a cloud with a linear 
size ~6pc and a mean density of 320 cm^'^. For these param- 
eters, the average column density is N{Il2)=5.9x IQ'^^ cm^^ 
and the corresponding average visual extinction is Ay =2.8 ™. 
The cloud is inhomogeneous and the extinction varies depend- 
ing on the line of sight. In the following, extinction and infrared 
spectra are only calculated toward one direction. The range of 
extinctions is Ay =0.31-24.6. Higher extinction means larger 
variations in intensity. In our library method, this should lead 
to larger errors, at least if the library always contains the same 
number of intensity intervals. 

We created two new models by scaling the densities of the 
original model. For the radiative transfer problem, the relevant 
parameter is the column density rather than the volume density 
and the results of the tests are equally apphcable to smaller and 



Table 1. Parameters of the model clouds used in the test- 
ing: name of the model, average density assuming linear size 
L=6 pc, and the average and maximum values of visual extinc- 
tion seen towards the direction for which the dust emission was 
computed. 



model name 



(n) (^v) Al 



[cm-^] 



CO 
CI 
C2 



320 2.76 
1280 11.04 
5120 44.16 



24.6 
98.5 
394.0 



correspondingly denser sources. The parameters of the models 
are listed in Table |l} The visual extinction is in all cases well 
below 1000™ and dust emission can be solved without iteration 
(see Sect. 2.1). 



3.2. The size of the library 

With discretization reduced to 32'^ cells, the emission spectra 
can be calculated directly within a reasonable time, and com- 
paris on w ith spectra computed with the method introduced in 
Sect. 2.2 makes it possible to directly estimate the eiTors caused 
by the library method. The radiation field is simulated using the 
same set of random numbers, and therefore the quoted errors 
reflect only the accuracy of the approximation. The incoming 
intensity is estimated with Monte Carlo methods, and the re- 
sulting random noise causes additional uncertainties. However, 
although the simulation of the radiation field represents a rather 
small part of the total computation time in our runs, the ex- 
pected noise at individual frequencies is only a few per cent. 
The effect on the dust emission is expected to be somewhat 
smaller, since the computation of dust temperatures involves 
averages over a number of frequencies. We computed maps of 
32 X 32 spectra using the 32^^ cell models and either direct so- 
lution or the library method. We first consider models without 
the PAH component, including only larger carboneous grains 
and silicate grains. 

With relatively low extinction, the model CO is the easi- 
est case. With two reference wavelengths, 0.12/im and 9.1/xm, 
both with 20 intensity intervals, the average relative error in 
the computed emission was below 2% in the range 5-300 fim. 
Taking into account all spectra and all wavelengths in this 
range, the maximum deviation between intensities computed 
directly and with the library method was ^6%. This is very en- 
couraging, considering that in the model the emission per unit 
volume (^density) varies over nearly three orders of magnitude 
and the intensity along different lines of sight over almost two 
orders of magnitude. The variations in the emission per dust 
grain are also significant, but much smaller, and this makes the 
use of the library method possible. 

The main computational cost lies in the solution of the dust 
temperature distributions, using the 20x20 spectra determined 
by the discretization of the intensity at the two reference wave- 
lengths. Even for such a small model, the computational time 
is significantly reduced, by a factor of ~ 30'^/20^ «68. 
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Fig. ^ shows the 100/im emission as the function of the lo- 
cal intensity at two reference frequencies for model CO. In the 
plot the contours are almost horizontal i.e. the emission can be 
predicted based on the intensity at only the shorter of the two 
reference frequencies. In models CI and C2 the higher extinc- 
tion makes the second reference wavelength more important. 
Fig. ^ shows clearly that in model C2 the spectrum of the local 
radiation field can no longer be approximated with only one pa- 
rameter (e.g. visual extinction), and the isocontours in the fig- 
ure are not parallel to either axis. At low intensities, i.e. deep 
inside the cloud, the emission follows the intensity at the longer 
wavelength and the contours are closer to vertical. At higher 
intensities (close to the cloud surface) the short wavelength ra- 
diation becomes relatively stronger and it is more important in 
determining the dust temperature and the infrared emission. As 
a result, the isocontours turn almost horizontal in Fig. ^. 

With the previous number of intensity intervals in the li- 
brary, 20^, the rms errors increase to 3.5% for model CI and 
6.2% for model C2. The errors might be reduced further by op- 
timizing the selection of the reference wavelengths. The wave- 
lengths should obviously be well separated from each other and 
coiTespond to sufficiently high optical depths so that intensity 
variations are well above the Monte Carlo noise. The optimal 
values depend on the radiation sources, the dust model and the 
cloud parameters. The reference wavelengths used above are 
based on some experimentation but are probably not optimal. 

In optically thick clouds, the inclusion of a third reference 
wavelength may be necessary. It is also better to use longer ref- 
erence wavelengths in order to avoid the extremely optically 
thick regime where accurate sampling of the radiation field 
becomes difficult. Using 20 intensity intervals at wavelengths 
0.20, 1.7 and 30yum, the rms errors are reduced to 1.8% and 
2.3% for models CO and CI, respectively. For the most opaque 
model, C2, the rms eiTor is 3.8%. It is, however, still possible 
to increase the size of the library by at least one order of mag- 
nitude, and the errors caused by the use of the library method 
can be kept below the noise due to the Monte Carlo sampling. 

The inclusion of PAHs does not cause significant problems 
for the library method. For the model CI, with a library with 
three reference wavelengths and 23^ intensity intervals, the rms 
error in intensities between 5 /im and 300 fim becomes 2.3%. 

It is not possible to give general estimates of the accuracy 
of the library method in more complicated situations. For ex- 
ample, the radiation field could probably not be estimated with 
sufficient accuracy using only two reference wavelengths if the 
model included two radiation sources with different emission 
spectra. The local field would clearly depend on the distance 
to both sources. Furthermore, in an inhomogeneous cloud the 
way each spectrum is attenuated depends on the density field 
which modifies both the intensity and the shape of the radiation 
spectrum. In principle, the accuracy can always be improved by 
increasing the number of intensity intervals and/or the number 
of reference wavelengths. The library method remains viable 
as long as the number of entries in the library is significantly 
smaller than the number of cells in the cloud. For example, 
with four reference wavelengths, each with 30 intensity inter- 
vals the number library entries is 8.1x10^ and method could 
still be used for all models with more that ~ 100'^ cells. 




-29 -28,5 -28 -27.5 -30.5 -30 -29.5 -29 



log(ii", refl) [cm Vn(H)] log(n'', refl) [cm Vn(H)] 

Fig. 2. The 100/im emission as the function of the local ra- 
diation field at reference wavelengths, Aiefi — 9.1/im and 
•^rcf2 = 0.12/im. Intensities are given as the number of pho- 
tons absorbed per cubic centimeter, divided by hydrogen den- 
sity. The upper frames correspond to model CO and the lower 
frames to model C2. 

3.3. Effect of cloud resolution on the created library 

When the library is first created, the model cloud is used only 
to estimate the range of intensities present in the cloud. In the 
library, the combinations of intensities at the reference wave- 
lengths do not correspond directly to any specific cells in the 
cloud. Therefore, the previous error estimates should also apply 
when the same library is used to calculate spectra from mod- 
els with higher discretization. The only exception may be the 
coldest spot of the cloud. Low resolution leads to smoothing of 
density peaks and as a result, a higher resolution model may 
contain a few cells where the intensity is lower than anywhere 
in the low resolution model. If the library does not contain en- 
tries corresponding to these lowest intensities the spectrum pre- 
dicted for the coldest core will have a spectrum that is too warm 
(i.e. the coldest spectrum of the low resolution model). 

There are several reasons why the effect is not expected to 
be very significant. The library represents a mapping between 
the intensity of the local radiation field and the dust emission 
per unit density, independent of the density. Comparing models 
with different discretization, the largest effect on the emitted 
intensity comes directly from the density. Furthermore, since 
the radiative transfer equation is integrated over the whole line 
of sight the relative statistical errors in the emerging intensity 
are smaller than the relative errors in local emissivity. 

We tested the effect explicitly with model CI discretized 
into 90'^ cells. The library was built as above (including PAH) 
but based on simulation of the radiation field in models with ei- 
ther 32^, 64'^ or 90'^ cells. Spectra computed with the different 
libraries were compared and the results are shown in Table. ||. 
The rms difference between the predicted IRAS intensity ra- 
tios, /eo/^ioo are below one percent and only in the coldest 
core at the position (-2.1, -1.1) the error increases to ^12%. 
The spectra for this position are shown in Fig. bI This shows 
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Table 2. Comparison of spectra calculated for the model CI. 
The cloud was discretized into 90"^ cells. The libraries were 
created using simulations of the radiation field in models with 
either 32'^, 64'^ or 90"^ cells. The last one was used also as the 
reference to which all other models were compared. The rms 
difference for range 10-300/xm as well as maximum deviation 
in IRAS ratio /go //loo are shown. The last column shows the 
comparison with direct calculations (see text). 



model 
rms(AI) 

max A(/6o//ioo) 



library (32'') 

0.4% 

12.4% 



lib. (64'') 

0.15% 

2.5% 



lib. (90'') 



direct 
1.3% 
1.4% 



^^ 

m 
m 




Fig. 3. Spectra calculated for the coldest core in model cloud 
CI using the library method. The libraries are based on simu- 
lations with models with either 32'^ (dashed line) or 90'^ (solid 
line) cells. This represents the overall worst case where the li- 
brary based on 32^^ cell model overestimates e.g. the IRAS ratio 
ho/hno by more than 10%. 

that the spectrum of the coldest spot in the cloud is not nec- 
essarily well reproduced when the library is based on a model 
with crude discretization. 

It is not possible to solve the dust temperature distribu- 
tions for all 90"^ cells in a reasonable time using direct cal- 
culations i.e. without the library method. The dust emission 
was solved, however, for two layers of cells (i.e. two times 
90 X 90 cells) that included the coldest core identified from the 
computed map of IRAS ratio /go / /loo- Altogether 2 x90 spec- 
tra were computed for different lines-of-sight along this shce. 
Comparison of this subset with spectra computed with the li- 
brary method shows rms difference of ^1-2%, irrespective of 
the library used. This indicates that, apart from the coldest spot, 
the libraries are equally good in predicting the infrared spec- 
trum. For the library based on simulations with a 90'^ cloud the 
maximum deviation at any single frequency is down to 5%. 

The effect of the cloud resolution on the created library is 
relatively small. There are at least two ways to reduce these 
errors further One can always create the library using the same 
discretization as used in the calculation of the spectra. This can 
be time consuming, but the library method will still be much 
faster than the direct solutions since the dust emission is solved 



Table 3. The range of surface brightness values (Jy sr^^) in 
IRAS bands for models discretized into either L'^=32'^ or 128^ 
cells. The power of ten is given in parenthesis. 



L 1(1 2^m) 



I(25^m) 



I(60Atm) 



I(lOO^m) 



32 
128 



6.6(4)- 1.1 (6) 
5.1(4)-2.0(6) 



1.5(5)-2.5(6) 
l.l(5)-3.5(6) 



4.3(5)-8.5(6) 
3.4(5)-1.6(7) 



4.8(6)-7.5(7) 
3.5(6)-1.4(8) 



only for the radiation fields selected in the library and not for all 
cells. Another possibility is to create the library using a model 
where the range of intensities is artificially increased. This can 
be achieved by creating the library using a model cloud with a 
slightly higher density. 

3.4. Effect of cloud resolution on spectra 

In this section we compare spectra calculated for model CI 
with cloud discretized into either 32'^ or 128^ cells. All calcu- 
lations are made using the same library which includes three 
reference wavelengths, each with 23 intensity intervals. With 
increasing resolution, the edges of high density regions are bet- 
ter resolved, and this might lead to changes in the intensity dis- 
tributions and/or intensity ratios. 

Table || shows the ranges of surface brightness values for 
the four IRAS bands for the two discretizations. As expected, 
with low resolution, the intensity variations are smoothed 
and the range of intensities is also correspondingly reduced. 
Observationally, this is more important at the limit of high sur- 
face brightness values, where the increase of the spatial resolu- 
tion by a factor of four increases the peak intensities by approx- 
imately the same amount. In the low-resolution model, each 
spectrum coiTesponds to an average surface brightness over a 
larger area, but the values themselves are not significantly bi- 
ased. For example, the predicted colour ratios shown in Fig. ^ 
are essentially the same irrespective of the resolution. The ratio 
Ii2 1 /25 changes slightly as the resolution is increased. This ra- 
tio traces the PAH component and higher values imply stronger 
radiation field. In particular, the values in the higher peak at 
I12/I25 ~0.55 originate mostly from the surface of the pro- 
jected cloud. 

The average properties of the emission, the surface bright- 
ness and colour ratios are almost independent of the discretiza- 
tion. Compared with the 32'^ model, however, the calculation 
with the higher resolution yields an approximately 2% lower 
total infrared luminosity. Although small, the change is sig- 
nificant, since the estimated uncertainty per spectrum is al- 
ready of the same order The most likely explanation is that 
at lower resolution, the holes that appear between the filaments 
are smoothed out, and this increases the total probability that 
incoming photons are absorbed within the cloud. Bernard et al. 



(1992) noted that for spherically symmetric models with power 



law density distributions, the approximation of the radial pro- 
file with a step function leads to Mnc/erestimation of the infrared 
emission. This effect must also apply to the emission predicted 
from the cores and filaments in three-dimensional models, es- 
pecially as the discretization is much more coarse by necessity. 
At least in the present case this did not seem to have affected 
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Fig. 4. Distributions of selected IRAS colour ratios in maps 
computed using either 32'^ (thin line) or 128'^ cell models (thick 
line). 



the average surface brightness values. However, without a rel- 
atively high resolution, the dense cores would not be resolved 
at all. 

Perhaps the main advantage of high resolution is that it en- 
ables statistical studies of the simulated infrared maps. In such 
studies the map size should preferably be at least 256^ in order 
to provide sufficient range of scales. 

4. MHD models 

In this section we present dust emission spectra computed with 
three-dimensional model CI consisting of 128'^ cells. Emission 
is calculated for t he ran ge 5/im-300//m using the dust model 
of Draine & Li (2001). The density field is the solution of 
the three-dimensional MHD equations, in a regime of forced 
highly supersonic and super- Alfvenic turbulence and with an 
isothermal equation of state (model C in Juvela et al. 2001). 



We first present results from the original model and compare 
them with calculations with lower discretization models. We 
also discuss the effect of the cloud structure on the observed 
surface brightness maps and, finally, study qualitatively the ef- 
fects introduced by spatially varying dust abundances. 

4. 1. The original model 

The library (i.e. the mapping between the incoming intensity at 
selected reference wavelengths and the resulting emission) is 
created using the original models resampled onto a 32^ cell 
grid. We use three reference wavelengths, each with 23 in- 



tensity intervals. As described above, the downsized model is 
needed only to estimate the range of intensities and profiles for 
the incoming radiation. Therefore, at this stage the spatial reso- 
lution is not crucial for the final accuracy. The use of the library 
leads to significant time saving since in the large, 128'^ cell 
model the radiative transfer problem is solved only at the three 
reference wavelengths. In this way the computational time is 
reduced by a factor of ^ 100 (ratio of total number of simulated 
wavelengths and number of reference wavelengths). The main 
cost lies, however, in the solving of the dust temperature distri- 
butions. Here the time saving for the 128"^ cloud is a factor of 
128'^/23'^ ^^170. As shown in Sect. ^ the expected rms error 
in the computed spectra is less than 2%. 

The surface brightness maps for the four IRAS bands are 
shown in Fig. || The maps have been calculated using the 
bandpass profiles given in the IRAS Explanatory Supplement 
(Beichmann et al. |1998 ). In all four maps the surface bright- 
ness probes the same filamentary density field. Maps of sur- 
face brightness ratios (Fig. ^ show, however, clear variations 
in the FIR colours and, as expected, dense cores and filaments 
are clearly visible through their lower colour temperatures. 

Figure ^ shows the average infrared spectrum and two indi- 
vidual spectra with extreme values of the ratio between 15/xm 
and 600/im intensities. The spectrum with relatively strongest 
mid-infrared emission originates from the upper left hand cor- 
ner i.e. gas subjected almost to the full intensity of the external 
radiation field. The coldest spectrum is observed at the core lo- 
cated at (-2.1, -1.1). There the emission peaks close to K:200/J,m 
but the core is already quite noticeable in the maps of IRAS 
colour ratios (Fig. ||). 

For the Li & Draine (2001 ) dust model, the ratio /eoZ-^ioo 



is close to 0.2 when dust is exposed to the pure interstellar ra- 
diation field. The dust model of Desert et al. (1990) predicts 
the same ratio Iqo/Iioo ^0.2. This is, of course, not surprising 
since both dust models are created to the reproduce existing 
IRAS and COBE measurements. The ratio of 0.2 is similar to 
those found towards many diffuse clouds (e.g. Laureijs et al. 



1991 ). This ratio applies, however, only to the diffuse medium. 
In our cloud CI, the average ratio is '^O. 1 1 . This cloud is quite 
opaque, with mean visual extinction of Ay ^11. The model 
CO has lower visual extinction, (Ay) ^2.8, and is closer to 
typical regions of molecular clouds. In CO, the average colour 
ratio is {leo/Iioa) "^0.15. This is the same as the value found 



by Abergel et el. (1994) for the Taurus molecular cloud com- 
plex. 

Some of the lowest colour temperatures are reached in the 
filament in the upper right hand corner of the maps. Like many 
other features, this is actually a superposition of separate fil- 
aments or sheet-like structures along the line of sight. In this 
direction, the densest filament is ~1 pc in diameter and has a 
maximum density of ^10^ cm~^. This roughly corresponds to 
a visual extinction of 10'" across the filament {^5™ to the cen- 
tre of the filament). There is, however, significant additional 
extinction from the more diffuse material between the dense 
regions. The ratio /eo/Iioo decreases within the filament by a 
factor of ~2.0 relative to the average value over the map, and 
is therefore approximately 4 times lower than the value given 
by the dust model in the diffuse medium. In observations the 
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ratio is often found to vary by a larger amount. In Lynds 134 
the ratio decreases by a factor of ~7 from the diffuse medium 
to the cloud cores, where the visual extinction is Ay ~10, 
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Fig. 6. Ratios between surface brightness values in IRAS 
bands. 
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close to that in our model CI (Laureijs et al. 1991). In the cirrus Z ' 



cloud MCLD 123.5+24.9 in the Polaris Flare the observed ratio 
/go/Zioo is 0.034 (Bernard et al. 1999) while the quoted aver- 
age v alue for cirrus clouds is (/eo/Iioo) ^0.31) (Lagache et al. 
1998 ). The ratio (/eo/Iioo) has therefore decreased by a factor 



of '^9, although the extinction is lower than in the filament in 
our model. As discussed by Bernard et al. ( 1999 ), the expla- 
nation of the observations requires spatial variation of the dust 
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properties. We will study these effects qualitatively in Sect. 4.3 



The calculations were caiTied out using a cloud divided into 
128"^ cells while the library was created with a model with only 
32"^ cells. However, according to Sect. ||the expected rms eiTor 
of the spectra is only ^2% and the error in /go /Iioa is less than 
15%. 

The average ratio (/i4o/^24o) calculated using the DIRBE 
system response functions is ~0.88 for model CI and ^1.05 
for mo del C O. In the spherically symmetric models of Bernard 
et al. ( 1992 ; see Lagache et al. 1998 Table 3), similar ratios 
are obtained only when the extinction is low compared with 
the average extinction in our models. For example, with Ay=4 
towards the cloud centre, their ratio has already decreased to 

0. 77, well below the value in model CI which has, never- 
theless, higher average visual extinction. This may be partly 
caused by differences in the dust models, but it is mainly the re- 
sult of the inhomogeneity of our models. In an inhomogeneous 
cloud, the ratio remains higher because there are no large, con- 
nected regions completely shielded from the external radiation 
field. 

The spatial power spectra for the 12/im and 100/^m maps 
are shown in Fig. || These are one-dimensional power spectra, 

1. e. they are computed from the coiTesponding two-dimensional 
power spectra by averaging over annuli of constant wavenum- 
ber. Although the MHD simulation used periodic boundary 
conditions, the external heating introduces radial gradients in 
the infrared maps, and therefore direct comparison with obser- 
vations is not possible at the largest scales. We will address this 




Fig. 7. Average infrared spectrum observed from model CI. 
The dotted lines show two individual spectra with extreme ra- 
tios between 15/im and 600/xm intensities. 



question in future work, where more detailed modelling will be 
done, also taking into account the asymmetry of the galactic ra- 
diation field. 

The finite spatial resolution of the MHD simulation and the 
relatively small size of the maps, 128 x 128 pixels, further limits 
the available dynamic range. In Fig.|]the vertical lines show the 
range where the spectrum follows an approximate power law, 
and where we determine the exponent of the power law. We 
obtain exponents /3 =-2.54 and -2.72 for 12/im and lOO/^m, 
respectively. In Sect. ^ we concluded the accuracy of the com- 
puted spectra is probably of the order of 2% per frequency 
point. Averaged over the IRAS filters the noise should be much 
lower. This is confirmed by the fact that the power spectra do 
not show a flat part typical of noise. Assigning an uncertainty 
of 2% to the computed IRAS values we estimated with Monte 
Carlo methods an uncertainty of ~0.01 for the exponent of the 
power law. This error estimate is valid for the comparison of 
maps with the same underlying density distribution. The expo- 
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nent could be changed significantly more by selecting a differ- 
ent range of spatial frequencies or e.g. a different projection of 
the cloud. 

The power law exponent for the column density map is - 
2.80. The same power spectrum slope is predicted analytically 
for the proje cted density in supersonic turbulence by Boldyrev 
et al. ( 2002 ), based on the v elocity structure function scaling 
predicted by Boldyrev ( 2002 ). 

Compared with the column density maps the infrared maps 
have more power on small scales but the slopes are in the range 
obtained from observations. The power law exponent of IRAS 
maps determined by Gautier et al. (1992) was k =-3.0. Based 
on ISO observations, however, Herbstmeier et al. (1998) found 
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10 



different power laws in different clouds, from very shallow 
spectra, fc ~ — 1, in dark regions up to k =-3.6 in the Cepheus 
flare. 

The power law exponent of the infrared maps is mostly a 
result of the underlying density distribution. The infrared maps 
do contain, however, significantly more power at high spatial 
frequencies, above the range used to compute the values of the 
power law exponents. This is due to dust temperature variations 
that cause intensity fluctuations on top of the column density 
distribution. The 12/im map seems to have more small scale 
structure than the 100/im one. At these high frequencies the 
underlying density distribution is affected by the finite reso- 
lution of the MHD simulations. However, for the given den- 
sity distribution we can compare the power spectra up to the 
Nyquist limit i.e. in Fig. ^ up to ^lOkpc"^. The differences 
were also observed in other models and they are well above 
the Poisson noise from the Monte Carlo simulation. Most of 
the fluctuations in the power spectra are caused by the underly- 
ing, common density distribution and e.g. the difference 12 /im 
and 100 /im is more significant than a comparison with the ap- 
parent noise in Fig. |] would lead to believe. The quantitative 
result could, however, change with better spatial resolution if 
the density field contained significant small scales structures. 

At 12/im, the emission mostly depends on the heating by 
short wavelength photons that are absorbed at the surface of 
the density structures. On the other hand, the heating of the 
large grains varies less rapidly, and one could expect the 100/im 
map to show relatively less small scale structures. Nevertheless, 
analysis of IRAS maps has not shown clear differences be- 
tween the different bands. It is conceivable that in the obser- 
vations, the differences are masked by projection effects or by 
the spatial variation in dust abundances and properties. 

4.2. Effects of inhomogeneity 

To check to what extent the three-dimensional cloud structure 
affects the properties of infrared emission, we compute spectra 
from two modified models, one with a smoother density struc- 
ture and one with a smaller density contrast than model CI. 

In the case of the first model, Cla, we re-arrange the cells 
in model CI so that in each column along the line of sight, the 
density increases monotonically towards a plane that is perpen- 
dicular to the line-of-sight and goes through the cloud centre. 
In this way, both the column density field and the probabil- 
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Fig. 8. Spatial power spectra of the 12/im (thin line) and 
100/im (thick line) surface brightness maps computed from 
model CI. The 100/im curve has been scaled down by a fac- 
tor of 5000. The Nyquist limit corresponds to fc = 10.7pc^^, 
but since the actual resolution of the MHD simulation is some- 
what lower, the power already drops above k ^^Spc^^. The 
dashed curve shows the power spectrum of the underlying col- 
umn density map with arbitrary scaling. The vertical lines limit 
the range used to compute the exponents of the corresponding 
power laws. 

ity density function of the volume density are the same as in 
the original model CI, but the density structure is smoother 
Comparison with the original model showed some noticeable 
changes in the surface brightness distribution, which now fol- 
lows more closely the column density distribution. In the orig- 
inal cloud, the observed infrared intensities are more sensitive 
to the intensity of the radiation field at the location of each 
individual filament or core. For example, the filament going 
through position (-2,2) is brighter in Fig. |] than in Fig. ^ be- 
cause it is located near the surface of the cloud. 

In model Cla, there are also some changes in the distri- 
bution of IRAS colours (see Fig. There is a tighter correla- 
tion between 60/im and 100/im intensities than in the original 
model CI. This is expected, as the cloud is inhomogeneous 
(clumpy or filamentary) essentially only in the plane of the sky, 
and smooth in the third spatial dimension. 

The modified cloud CI is not 'clumpy' along the line 
of sight, and following the preceding discussion on the ef- 
fects of resolution, one might expect higher infrared emission. 
However, the emission is some 20% lower. This is a direct 
consequence of the fact that the modified cloud is effectively 
smaller along the line of sight. Since most of the mass is con- 
centrated to a smaller volume, less radiation enters dense parts 
of the cloud and the total absorbed energy is also coiTespond- 
ingly lower Conversely, if dust properties and the external radi- 
ation field were known to an accuracy of a few percent, it would 
be possible to make rough estimates of the line-of-sight size of 
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Fig. 9. IRAS colours for model CI a, a modification of model 
CI where the density structure along the line of sight has been 
made smooth by re-arranging the cells along the line of sight 
(see text). 

clouds, even when the degree of dumpiness in this direction is 
not known. 

In model CI, the ratio between the highest and lowest den- 
sities is close to 6x10"*. We produced another, more uniform 
cloud model, Clb. By applying a linear scaling to the density 
values, the total density contrast was reduced by a factor of 100. 
The average column density remains the same as in models CI 
and CI a, but the absorbed energy is higher since it depends on 
column density via an exponential function. As a result, com- 
pared with the original model CI, the total emitted energy in- 
creases by ^7%. For the IRAS bands, the difference is slightly 
higher, ^^8.5%. In spite of the increased energy input, the av- 
erage ratio between IRAS bands gives a lower colour tempera- 
ture. The difference is probably real, although it is only a cou- 
ple of per cent, which is close to the expected accuracy of the 
calculations. In model Clb, near- and mid-infrared emission 
is concentrated on the surface of the model cloud, and the in- 
ner parts are correspondingly cooler The infrared colours were 
averaged over all spectra without weighting them with the cor- 
responding intensities. Consequently, the average colours can 
indicate lower dust temperatures than in the model CI, where 
the short wavelength emission is spread more evenly over the 
projected cloud area. 

The power spectra of the 12/im and 100/im maps are al- 
most identical in the range (see Fig. ^ used for the determi- 
nation of the power law exponent and the exponent is ^-2J 
for both maps. If the power spectra are scaled on top of each 
other in this range of spatial frequencies, the 12/im curve is 
only sUghtly above the 100/im curve at high frequencies, while 
at low frequencies, the 12/im is significantly higher, due to the 
stronger brightening towards the cloud surface. 

4.3. Dust abundance variations 

There is strong observational evidence that the dust properties 



change, not only in dense cloud cores (e.g. Laureijs et al. 1996) 



but also in clouds with moderate visual extinction (Bernard et 



al. |1999| ; Cambr esy e t al. |2004 Stepnik et al. |2004 Miville- 
Deschenes et al. 2002). Possible causes are the growth of grain 
mantles and the formation of large dust aggregates. The pro- 
cesses can take place only in dense medium that is protected 
from the external UV radiation field. 

We use the three-dimensional model CI to study the qual- 
itative effects of abundance changes on the observed infrared 
radiation. The model has mean density ~ lO'^ cm^'^ with peak 
values higher by two orders of magnitude. 

We adjust the abundances of the grain populations accord- 
ing to the local density only. The modified model is called Clc. 
At each position, the dust is composed of five components, 
the first three of which ar e PAH , silicate and carbon grains 
of the Li & Draine model (2001). The additional two compo- 
nents are modified silicate and carbon grain populations with 
a flatter size distribution: the original exponent of the size dis- 
tribution is adjusted so that the original number ratio of grains 
below and above 25 nm is reduced by 50%. At low densities, 
the abundance of the new components is zero. Their abun- 
dance increases linearly with the logarithm of the density. The 
abundance of the unmodified grain populations is decreased 
linearly with respect to the logarithm of density and the total 
mass of dust remains unchanged. Above 10^ cm^"^, all dust is 
in the modified silicate and carbon grain populations. This is, 
of course, only a qualitative description of the expected dust 
abundance variations. In particular, we ignore all changes in the 
dust optical properties expected to be caused by the formation 
of aggregate grains. As a result, the model can be expected to 
underestimate true colour temperature variations across dense 
filaments or cores (see e.g. Stepnik et al. 2001 ). 

Figure ^ shows maps of IRAS intensity ratios in model 
Clc. Since we have reduced the number of small grains the av- 
erage colour temperatures have dropped. The ratio {I25/I10Q) 
has decreased the most, by 34%. The ratio {Iqq/Iioo) is now 
lower by 23%. The average 100/im surface brightness has only 
increased by and the 60/im emission is below that of 

model CI. 

The 12/im map contains relatively more power at small 
scales. The short wavelength maps are not expected to change 
much, since most of the emission comes from low density 
material where grain abundances remain almost unchanged. 
Abundance variations should increase far-infrared emission in 
the densest regions, and the corresponding maps are, therefore, 
expected to contain more power at high spatial frequencies. 
However, the changes in the power spectra remain very small. 
The new exponents of the power laws are -2.51 and -2.67 for 
12/im and 100/im, respectively; the 100/im slope has only flat- 
tened by '--^2% relative to the original model CI. 

In the filament at position (2,2) the 12/i,m surface bright- 
ness has decreased by ~30%, while the 100/xm emission has 
increased less than 10%: this is less than the average increase 
over the whole map. The change in the dust abundances leads to 
lower extinction, especially at short wavelengths, and this de- 
creases temperature variations across the filament. The lowest 
values of the ratio /go/^ioo are ^0.05 in this region, i.e. some 
30% lower than in the original cloud CI. Compared with the 
diffuse medium value of 0.2, the reduction in /go/^ioo is still 
only a factor ~4. This is still less than found in many observa- 
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Fig. 10. Maps of selected IRAS colours calculated for the 
model Clc. In the model the dust abundances vary depending 
on the local density (see text). The colour scales are the same 
as in Fig.0. 



tions (e.g. Laureijs et al. 1991 ; Bernard et al. 1999). The model 



cloud has very significant extinction, even outside the densest 
filaments, and this should help to further decrease the /go / /loo 
ratio. 

In the presented model, we have only adjusted the grain 
size distributions by halving the ratio of grains below and 
above 25 nm in dense cloud regions with n > 10"' cm^'^. There 
is some evidence that dust properties can undergo significant 
changes even at much lower densities, below 10"^ cm^"^ (e.g. 



Bernard et al. 1999). On the other hand, sub-mm observations 



of dense cores seem to indicate that in those environments, all 



small grains may be removed (e.g. Stepnik et al. 2001 ). If small 



dust grains stick onto the surfaces of larger grains, the optical 
properties of the grains will change and lead to a lower phys- 
ical temperature. As a result, the far-infrared colour temper- 
atures decrease, and this can explain the very low /eo/Aoo- 
ratios observed. In dense environments, dust may also form 
larger, fluffy aggregate particles that have similarly enhanced 
far-infra red em issivity and low p hysical temperature (Mathis & 
Whiffen |l989t Ossenkopf [l993t Wright |1987|; Bazell & Dwek 



199C). Many of the observed cold cores are too small for the 



low colour temperatures to be explained merely by the attenua- 
tion of the external radiation field. The need for large variations 
in the abundance ratio of large and small dust grains is, there- 
fore, already well established. 

5. Conclusions 

We have presented a new approximate method that can be 
used to compute the infrared dust emission from large, three 
dimensional cloud models consisting of millions of compu- 
tational cells. The method is based on simple discretization 
of the incoming intensity using a few reference wavelengths. 
We have demonstrated that relative accuracy of a few per- 
cent is easily reached, even when transiently heated grains are 
included. Such accuracy is hard to match either by observa- 



tions or the Monte Carlo sampling that is often used to es- 
timate the intensity of the radiation field. We have tested the 
method using model clouds with average visual extinctions up 
to Ay '-^lOO. However, with use of more reference wavelengths 
and/or higher discretization, the method can be applied to even 
denser clouds, and to models with several heating sources with 
different emission spectra. 

We have used our new approximate method to compute 
dust emission, from near-infrared to sub-mm wavelengths, 
from models based on MHD simulations with supersonic and 
super- Alfvenic turbulence. The spatial discretization was found 
to have little impact on the derived distributions of the intensity 
ratios between IRAS bands. On the other hand, spatial averag- 
ing will always reduce the range of intensities. High resolution 
calculations are therefore necessary, especially when the com- 
puted infrared maps are used for statistical studies or for direct 
comparison with observations. 

The spatial power spectra at different wavelengths gener- 
ally show little variation. The shorter wavelengths have, how- 
ever, slightly higher relative power at higher spatial frequen- 
cies. For the selected MHD model, the power spectra follow 
approximately a power law with an exponent between -2.5 and 
-2.7, while the power spectrum of the underlying column den- 
sity map has a slope of -2.8. 

We qualitatively studied the effect of density dependent 
dust size distributions on the observable infrared maps. By re- 
ducing the relative abundance of very small grains by 50% in 
regions with densities n > lO^cm"^, the ratio len/Iioo could 
be brought down to ^0.05 in dense regions where visual ex- 
tinction is Ay ~10. In observations, the colour ratio can reach 
lower values even in less opaque clouds. This indicates that 
dust properties must undergo even larger changes, with the re- 
moval of most of the small grains and the growth of larger ones. 

The numerical method presented in this work will be very 
useful to constrain dust properties by comparison with new ob- 
servational data. Dust column density maps of large molecular 
cloud complexes can now be obtained from stellar extinction 
measurements based on the " Two Microns All Sky Survey" 
(2MASS). Padoan et al. (2002) have recently generated extinc- 
tion maps of the whole Taurus molecular clouds complex us- 
ing 2MASS data, with a spatial resolution comparable or better 
than that of 100 /im IRAS images, spanning a range of visual 
extinction from Ay « 0.3 mag to Ay « 30 mag. The com- 
parison of this type of dust column density maps with IRAS 
images and, in the near future, with SIRTF-MIPS images will 
provide new observational constraints over a large range of col- 
umn density and gas density values. Detailed three dimensional 
calculations of dust emission from physical models of interstel- 
lar clouds will be necessary for reproducing the observational 
data and improve our knowledge of the properties of interstellar 
dust. 
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